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Abstract 

Nonthermal scaling phenomena can exhibit a characteristic 
dependence on the dimensionality d of space. For d = 3 and 4 we 
simulate a relativistic scalar field theory on a lattice and compute 
turbulent scaling exponents. We recover Kolmogorov or weak wave- 
turbulence in the perturbative high-momentum regime, where it 
exhibits the scaling exponent = d — 3/2. In the nonperturbative 
infrared regime, we find a different scaling exponent Kg = 4 (5) 
for d = 3(4), which is in agreement with the recently predicted 
anomalously large values = d + 1 oi strong turbulence. We show 
how the latter can be seen to characterize stationary transport of a 
conserved effective particle number. 



1 Introduction 



Turbulent scaling phenomena appear for very different physical systems 
and length scales, ranging from early universe reheating dynamics [H [2] to 
physics of supernova explosions [3] or laboratory experiments with ultra- 
cold quantum gases [5]. While many aspects of turbulence have long 
reached textbook level rather little is known about turbulent behavior 
in nonperturbative regimes of relativistic quantum field theories. Here the 
strong interest is to a large extent also driven by related questions concerning 
relativistic heavy-ion collisions |6]. 

It has recently been demonstrated that a new class of turbulent scaling 
phenomena exist in the nonperturbative regime of relativistic scalar field 
theories for sufficiently high occupation numbers per mode |2]. The 
nonthermal scaling solutions were shown to exhibit a strong enhancement 
~ |p|~^ in the infrared for three spatial dimensions, while at high momenta 
a well-known scaling regime of weak wave-turbulence ~ |p|^^^^ is observed. 
In Ref. [7j a characteristic dependence of the nonperturbative infrared 
solutions on the dimensionality of space was predicted. Turbulent scaling 
phenomena are insensitive to the details of the underlying microscopic 
theory. In particular, the infrared scaling solutions are expected to belong to 
universality classes. These may only depend on few general properties such 
as space dimension, symmetry, field content and conserved charges. 

In this work we compute turbulent scaling behavior of a self-interacting 
iV-component scalar field theory. This is done using classical-statistical 
simulations on a lattice, which are expected to be an accurate description also 
for the corresponding quantum theory for the considered high occupation 
numbers per mode [2]. In order to further classify the universality class 
of nonperturbative scaling solutions, we perform simulations in three and 
four spatial dimensions. In addition, we simulate for = 4 and N = 10 
field components. We find striking agreement of the numerical results 
with the previously predicted analytic dependence on the dimensionality of 
space based on resummed large- techniques |7]. In particular, we find 
no indication for a dependence on and the results are consistent with 
a vanishing anomalous dimension as well as a dynamic scaling exponent 
z = 1 for the relativistic theory. We show that the phenomenon of strong 
turbulence at low momenta may be associated to stationary transport of a 
conserved effective particle number. 

The paper is organized as follows. We start with the definition of suitable 
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correlation functions, which can be used to discuss weak as well as strong 
turbulence from a common framework in Sec. [21 After a brief review of 
perturbative wave-turbulence using kinetic theory, we analytically discuss 
strong turbulence as the stationary transport of a conserved effective particle 
number. In Sec. [3] we discuss the numerical results from lattice simulations. 



2 Stong versus weak stationary turbulence 

We describe turbulence in a relativistic self-interacting scalar field theory 
as the stationary transport of a conserved quantity. Typically this is 
done perturbatively using kinetic theory. To obtain a description that is 
valid also beyond perturbation theory, we first define correlations functions 
which suitably describe the physics of turbulence. These are evaluated 
perturbatively using kinetic theory in Sec. 12. II and then nonperturbatively in 
Sees. O and El 

For a real scalar field theory with Heisenberg field operator 
where x = x) denotes the time x° and (i- dimensional space variable 
X, we consider the expectation value of the anti-commutator F{x,y) and 
commutator p(x, y) of two fields, 

F(x,y) = i({$(x),$(y)}) , p{x,y)=tmx),<l>{y)]), (1) 

respectively. Here the real spectral function p{x, y) is related to the retarded 
or advanced propagators by 

y) = G\y, x) = e(x° - p{x, y) , (2) 

with the step function Q{x^ — y^) = 1 for x° > and zero otherwise. 
Because of anti-symmetry, p{x,y) = —p{y,x), one has p{x,x) = 0. The 
real statistical two-point function is symmetric, F{x,y) = F{y,x). It is 
convenient to introduce Wigner coordinates 

X^ = ^^ , s^ = x^-y^ (3) 

and Fourier transform with respect to the relative coordinates, 

Fp{X) = j d-^+^s e-'P^'"F{X + s/2, X - s/2) , (4) 

p,(X,^-./d«.e--V(X + .V2.X-./2) (5) 
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with p = {p^,p). The conventional factor of —i in the transform for the 
spectral function makes the latter real in Fourier space and we use a tilde to 
denote this. According to ([2]), retarded and advanced propagators in Fourier 
space are then related to the spectral function by 

tp,{X) = G^{X) - G^{X) . (6) 

For the real scalar field theory one finds from the definitions ([1]), 

= F,iX) , p.,{X) = -p,{X) . (7) 

Without loss of generality we write 



Fp{X) = (n,(X) + ~p,{X) 



(8) 



which defines the function np{X) for any given Fp{X) and pp{X). We 
emphasize that without additional assumptions ([8]) does not represent a 
fluctuation-dissipation relation, which holds only if np{X) is replaced by a 
thermal distribution function. We will not assume this in the following and 
keep np{X) general at this stage. In particular, (JTj) then implies the identity 

n.p{X) = - (ripiX) + 1) . (9) 

Since we are interested in stationary behavior, we may consider the 
dynamics employing a gradient expansion to lowest order in the number 
of derivatives with respect to the center coordinates X'^ and powers of the 
relative coordinates s^. This is a standard procedure for the derivation of 
kinetic equations from field theory, and for the spectral function to lowest 
order one has [HI E] 

'2jf^^Pp{X) = . (10) 

For spatially homogeneous ensembles (fTO|) implies a constant pp that does 
not depend on time. In contrast, the statistical function Fp{t) to this order 
can depend on time t = X° and using ([8]) we can write 

r r dp' , dupit) _ 

i W/'^'^ - j, ^^PPv^^ -C[n]{t^v)- (11) 

Here C[n] denotes the 'gain' and 'loss' terms, which describe the effects of 
interactions to lowest order in the gradient expansion. We will determine 
C[n] using perturbation theory in Sec. 12. ![ and in Sec. I2.2l it will be obtained 
from a resummed expansion to next-to-leading order (NLO). 
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2.1 Weak wave-turbulence 



In this section we review some relevant aspects of perturbative Kolmogorov 
or weak wave-turbulence [10], which will be used below for comparison with 
the nonperturbative regime of strong turbulence. The free spectral function 
is 

pI = 2nsgnip^)5{{py-col) (12) 

for a relativistic scalar field theory with particle energy Up. Plugging this 
into f lTT]) gives 

^ = C[n]it;p), (13) 

where 

upit)^ I -£2p'pln,{t). (14) 

For a compact notation we will frequently suppress the t-dependence and 
write Up. In kinetic theory, when two particles scatter into two particles, the 
collision integral on the RHS of f|T3|) is of the form 

C^^\p) = J dn^^\p, 1, q, r) [(1 + np){l + n,)n^nr - npn,{l + n^){l + rir)] 

(15) 

where the gain and loss terms in the integrand take into account that for 
bosons there is an enhancement of the rate if the final state is already 
occupied. The details of the model enter J dil^^^(p, 1, q, r), which for an 
0(A^)-symmetric scalar field with quartic A/(4!A^)-interaction reads [H] 

I dn'^'ip, 1, q, r) = ^ (27r)'^+i 5^''\p + I - q - r) 

with the notation = J d'^p / [2tt)'^ . Using the approximation f[T^ with 
f|T6|) for the collision integral on the RHS of f|T3|) . one obtains the well- 
known Boltzmann equation for a gas of relativistic particles. Clearly, this 
approximation cannot be used if the occupation numbers per mode become 
large and, parametrically, for a weak coupling A a necessary condition is 
Up -C 1/A [TT]. In Sec. 12.21 we will discuss suitable approximations that are 
valid also for nonperturbatively large occupation numbers. 



4 



Turbulence is expected for not too small occupation numbers per mode 
such that quantum corrections can be neglected. For a regime 1 ^ rip ^ 1/A 
one may still use the above Boltzmann equation, which becomes 



dt 



J dri^^^(p, I, q, r) [(np + ni)nqnr - raprai(raq + Ur)] 



(17) 



For the considered real scalar field theory the energy density e is conserved. 
Since we restrict our discussion in this section to number conserving 2 o 2 
scatterings also the total particle number ntot is conserved, and 



Wp^p , ntot = Up. (18) 
Jp 

The fact that they are conserved and interactions are sufficiently local in 
momentum spacqj may be described by a continuity equation in Fourier 
space, such as 

^(c^pnp) + Vp-jp = (19) 

for energy conservation [10]. Similarly, particle number conservation is 
described by formally replacing Wp — > 1 in the above equation and 
corresponding substitution of the fiux density. For isotropic ensembles we 
consider the energy fiux A{k) through a momentum sphere of radius k. Then 
only the radial component of the fiux density jp is nonvanishing and 

/'Vp-jp = [ jp-dAp = {2nfA{k). (20) 

J p J dk 

Since Up is constant in time, we can thus write with the help of (HM 

For stationary turbulence the fiux A{k) is scale independent, i.e. the 
respective integral does not depend on the integration limit k. We consider 
scaling solutions 

Up = s'^'^Usp , Up = s^'^Usp, (22) 



^Strongly non-local contributions are suppressed by phase space and energy- momentum 
conservation. The classical theory requires, of course, an ultraviolet cutoff to regularize 
the Rayleigh- Jeans divergence. 
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with occupation number exponent k„ and assuming a linear dispersion 
relation. Since the physics is scale invariant, we can choose s = l/|p| such 
that Hp = IpI"'^"' rii and Up = |p| wi- Using the scaling properties f l22|) one 
obtains from f lT6|) 

J dn2^2(p,l,q,r) = s"^* J dn^^\sp,sl,sq,sr). (23) 

Here the scaling exponent fi4 for the theory with quartic self-interaction is 
given by 

= (^3d - A) - {d + 1) = 2d-5, (24) 

where the first term in brackets comes from the scaling of the measure and 
the second from energy-momentum conservation for two-to-two scattering 
described by (fTSjl . Apart from the quartic self-interaction, it will be relevant 
to consider also scattering in the presence of a non- vanishing field expectation 
value such that an effective 3- vertex appears pTj. In this case, one obtains 
along these lines 

/i3 = (2rf-3) - (rf + 1) = (25) 

To keep the discussion more general, we write for the scaling properties of 
the flux for a given m-vertex 

A{k) = --——1- rdp|p|^-i+i-''»('"-^)+'^-a;i^. (26) 



2^7r<^/2r(d/2 + 1) y ' dt 

For m = 4 this can be directly verified to agree to the two-to-two scattering 
case with a 4- vertex using (ITTl) and (IMl) . If the exponent in the integrand of 
f l26|) is nonvanishing, the integral gives 

Mk) -r4 7 71 ^1^- (27) 

Scale invariance up to logarithmic corrections is, therefore, obtained for 

d+1- K„{m-1)+ fXrn = 0- (28) 

This yields for the energy cascade 

m=4 , 4 m=3 , 3 

= d — - , K„ = a— -. (29) 
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One observes that stationary turbulence requires the existence of the hmit 



drill dt 



const 7^ 



(30) 



(i+l-Kw(m-l)+Atm-^0 d+1 — K„(m — 1) + /i. 



such that the colhsion integral must have a corresponding zero of first degree. 
Similarly, starting from the continuity equation for particle number one can 
study stationary turbulence associated to particle number conservation. For 
instance, for two-to-two scattering this leads to = — 5/3, and = d — 2 
for the case of interaction through a 3-vertex. 

2.2 Strong turbulence 

The above perturbative description of stationary turbulence becomes invalid 
at low momenta |p|, since the occupation numbers Up ~ Ipl"**" can grow 
nonperturbatively large in the infrared. This concerns positive values of 
the scaling exponent given by ( l29l) . which is the case for the dimensions 
d = 3 and 4 to be considered below. To understand where the picture of 
weak wave-turbulence breaks down and to compute the properties of the 
infrared regime, we have to consider nonperturbative approximations. In 
this section we consider the expansion of the two-particle irreducible (2PI) 
effective action in the number of field components to NLO to get analytical 
insight [12]. We extend the discussions of Refs. [21 E] by showing that a 
conserved effective particle number characterizes strongly modified scaling 
properties in the nonperturbative low-momentum regime. 

At NLO in the 2PI expansion the evolution equation (ITT]) for 




(31) 



reads 




C^^""N(t;p). 



(32) 



at 
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Here the NLO contribution [T2] can be written as 



CNLO(p) = j an^^^p, I, q, r) [(1 + np)(l + n^n^n, - npn^l + ng){l + n,)] 

+ J dfi(^^(p, /, q, r) [(1 + np)(l + ni)(l + nq)n^ - npninq{l + n^)] 

+ j dnl^^{p,l,q,r) [{1 + np)ningnr - np{l + ni){l + ng){l + Ur)] 

+ y dJl°'^^(p, /,g,r) [(1 + np)(l + ni){l + ng){l + rir) - npniUqUr] 

(33) 

where we consider the case of a vanishing macroscopic field, i.e. ($) = 0, 
relevant for the infrared regime [2]. Again we suppress in the notation the 
time dependence. Here 

/■jo2«2/ , N ^ p d/d/OdgOdrO [ , 

X pppipgpr [Ksip + + ^csip -q) + Ksip - r)] , 

/■ ,r.io3. , ^ ^ p d/d/OdgOdrO /■ ,^,1^, 

X PpPiPqPr [Aefr(p + + -^ofr(p + + Aefr(p - r)] , 

X PpPlPqPr KsiP - I) , 

/'^oo«4/ , N ^ p d/d/°dgOdrO r , 

X ppplpqpr Ke{P + I) ■ (34) 

In contrast to (ITB]) . the above expressions still contain the integrations over 
frequencies and spectral functions. We emphasize that the latter are, in 
general, not of the free field form f|T2|) in the nonperturbative regime. No 
quasi-particle assumptions has been employed and the only approximations 
are the expansion to NLO and the gradient expansion underlying ( ITTi) . 

The effective momentum- dependent 'coupling' Aefr(p) appearing at NLO 
in the 2PI expansion is given by [21 [7] 

W) = ^^j^, (35) 
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which involves the squared absolute value ll+n-^d))!^ = [l+n^(p)][H-n'^(p)] 
of the retarded or advanced self-energy 



R,A 

q 



(36) 



One observes that for sufficiently large p, for which 11^'^ (p) ^ 1, the effective 
coupling fl35|) approaches A. In this case the '2 <H- 2' contribution of the first 
line in ( 133|) is reminiscent of the two-to-two scattering process described 
by the perturbative expression presented in f lTSj) . The main difference is 
that the latter assumes a S-\ike spectral function such that all momenta are 
on shell. Therefore, in the perturbative expression f[T5]) off-shell processes 
involving the decay of one into three particles or corresponding 3 — )■ 1 
annihilation processes or even 0^4 processes are absent. They can occur 
in principle at NLO in the 2PI expansion, which leads to the different 
terms contributing to the RHS of fl33p . At sufficiently high momenta these 
off-shell contributions should be suppressed along with quantum-statistical 
corrections such that the spectral function approaches a (5-like behavior. In 
this case we would recover stationary turbulence characterized by a weak 
scaling exponent described above. In contrast, in the infrared Xcs{p) may 
have a nontrivial momentum dependence, which is discussed in the following. 
Following similar lines as in Sec. 12. 1^ we look for scaling solutions where 
1 such that (!33|) can be approximated by its classical-statistical limit 

Cd^^(p) = j dn'^^^{p,l,q,r)[{np + ni)nqnr - npni{nq + Ur)] 
+ / dnl^^{p, I, g, r) [{up + ni)nqnr - npni{nq - n^)] 

. jnM(n,^n.)n,n^,nMn,^n^)]. (37) 

Again, only the first term above is reminiscent of the perturbative expression 
in f lT7|) . In principle, nonperturbative scaling phenomena may involve an 
anomalous scaling exponent for pp = p{p^,p). Using isotropy we write 
following Ref. [7] 

p(/,p) = s2-'^p(,y,sp), (38) 
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with a nonequilibrium 'anomalous dimension' rj. A dynamical scaling 
exponent z is taken into account since only spatial momenta are related 
by rotational symmetry and frequencies may scale differently because of the 
presence of (non-)thermal corrections. Scaling behavior of the statistical 
correlation function 

F(/,p) = 5^+"^ F(sy , sp) (39) 
then translates with ([8]) for Up':^ 1 into 

n{p^,p) = s^'^+^nis'p^sp) . (40) 
With this one can determine the scaling behavior of Xcs{p)- From fl36|) follows 

n^'^(/,p) = s^n^'^(sy,sp) (41) 

with 

A = 4-d-z + Ks-'n- (42) 
If A > one finds from (155]) the infrared scaling behavior 

Aeff(/,p) = s-2AA,g(sy,sp). (43) 

For A < the effective coupling becomes trivial with Xes{p) — A, on which 
we comment below. Using these scaling properties one obtains from (|3] 



(44) 

and the same scaling behavior for J dQ^^^{p,l,q,r), J dr2^^'^(p, /, g, r) and 
JdQ^^\p,l,q,r). 

Similar to Sec. 12.1] for any conserved quantity we can compute the flux 
through a momentum sphere k. Stationary turbulence solutions then require 
that the respective integral does not depend on k. Obviously, energy is 
conserved. The highly nontrivial question is whether a conserved effective 
particle number exists since there is no conserved charge associated to particle 
number in the real scalar field theory. Off-shell processes included in (!37|l . 
such as 1 ^ 3 processes, make this manifest. In the following, we analyze 
whether the effective particle number nefj(t;p) given by (1311) represents a 
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conserved quantity for the nonperturbative low-momentum regime. Similar 
to fl21l) . the flux for this effective particle number now reads 



k 



dncs{t; p) 



<=ttV ; 2'^7r'^/2r(rf/2 + 1) y dt ^ ' 

The momentum integral can be evaluated along the lines of Sec. 12.11 using 
the above scaling properties with 

dt 'P' dt ' ^ ^ 

such that 

if the exponent in the integrand is nonvanishing. Scale invariance up to 
logarithmic corrections may, therefore, be obtained in the nonperturbative 
low-momentum regime for 

Ks = d + z -7]. (48) 

This scaling solution is associated to a conserved nes{t', p) for sufficiently low 
momentum pH Similarly, for the scaling solution associated to conserved 
energy one finds, taking into account an additional power of p^, the exponent 
Kg = d + 2z — 7] in accordance with Ref . [7] . 

The above discussion shows that in the presence of a conserved ngg (t; p) 
there is a strongly modified infrared scaling behavior as compared to 
perturbative treatments. In particular, fj48|) predicts a characteristic 
dependence on the dimensionality of space d and no dependence on the 
number N of field components. 



3 Lattice simulations 

In this section we solve the evolution equations for our theory in the classical- 
statistical limit using simulations on a lattice. Varying the dimensionality of 

^Refs. [ZIH assume that the 1 O 3 and -H> 4 contributions in ([57]) vanish to obtain the 
solution (|48| . See also the discussion of this point in Ref. [5 . We note that the momentum 
integral over neff(p) can be strongly infrared divergent for the discussed scaling solutions 
and requires an infrared cutoff. 
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space d and the number of field components N, we will then compare the 
numerical results with the analytical estimates f l48p in the infrared and f l29p 
for high momenta. 

The field equation of motion for the classical A^-component scalar field 
theory reads with 1 < a < A^: 

A ^ 

(^a{t, X) = (Ax - m^) Lfait, x) - — ^ ^bit, x)ipb{t, x)v9a(t, x) . (49) 

6=1 

For the numerical implementation of the above equation the leap-frog 
discretization is used on a cubic space-time lattice in three and four spatial 
dimensions. The initial conditions are chosen such that the system will 
evolve closely to non-thermal scaling solutions. To achieve this one can start 
with a nonequilibrium instability, such that low-momentum modes get highly 
populated [2]. Such instabilities are, for example, the tachyonic instability or 
the parametric resonance instability, which also have cosmological relevance 
as models for reheating [T3l [TH [15]. In this study we use initial conditions 
triggering parametric resonance: the space average of the field has a nonzero 
initial value {(pi{t = 0)) = 0o while {(pait = 0)) = for 1 < a < N. 
The nonzero momentum modes are initialized with a small amplitude white 
nois^ to provide a seed for unstable modes. The results are then averaged 
over different realizations of the initial noise distribution. 

Our main observable, the momentum dependent particle number is 
defined by 

1 ^ 

a=l 

where (pa{t,p) is the spatial Fourier transform of the field in d dimensions, 
i.e. (pa{t,p) = 1/VV J d'^xv9a(t, x) exp(ipx) with the spatial volume V. 

Fig. [T] shows the particle number spectrum for a three dimensional 
simulation using a 192^ lattice with A = 24, = and = 4. The infrared 
modes exhibit a slow time evolution, whereas higher-momentum modes seem 
to settle much more quickly. For the final plotted time (pot = 184000 one 
observes two separate regions with clear power laws. For high momenta 
the scaling exponent ^ 1.5 agrees well with the analytic prediction 

■^Its spectral composition is not important as long as the amplitude is small. 
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Figure 1: The particle number spectrum of the d = 3 simulation for different 
times in units of the initial field amplitude (pQ. 
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for Kolmogorov turbulence in three space dimensions as expected [TU]. 
As the occupation number per mode grows towards lower momenta, the 
perturbative approximation breaks down for the description of the infrared 
modes. Accordingly, one observes a strongly modified power-law. Assuming 
z = 1 for the relativistic theory and rj = 0, which is the case also at high 
momenta, the observed value ~ 4 agrees well with the estimate fllS]) in 
accordance with the results of Ref . [2] . 

A crucial test for the interpretation of these results in terms of the analytic 
estimates of Sees. 12.11 and 12.21 is their predicted characteristic momentum 
dependence. In Fig. |2]we show the particle number spectrum for simulations 
in d = 4 for different initial conditions and coupling values. The upper graph 
indeed shows a low-momentum scaling exponent ~ 5 as well as a high- 
momentum scaling exponent k^^ ~ 2.5, which are in remarkable agreement 
with the predicted values fj48|) for z = 1 and 1] = as well as f l29|) . The lower 
graph shows results for a much higher energy density in lattice units OLat- 
While the low-momentum scaling behavior is insensitive to these changes, 
only the upper graph shows the perturbative Kolmogorov scaling at high 
momenta. In particular, the observed high-momentum behavior is closer to 
the classical thermal exponent value of one, rather than to the Kolmogorov 
exponent for d = 4. 

We have seen that, in particular, the strong turbulence regime is very 
insensitive to details of the underlying theory such as couplings or initial 
conditions. It remains to see whether there is a dependence of the scaling 
behavior on the number of field components A^. In Fig. [3] we show results from 
simulations using A^ = 10 fields in three dimensions. Again the exponents 
follow very closely the analytic estimates for z = 1 and rj = 0. This 
indicates that the universality class for the turbulent scaling exponents does 
not depend on A^ in accordance with the analytic estimates. 

4 Conclusions 

Stationary turbulence is associated to conserved quantities. We have demon- 
strated that the nonperturbative scaling solution first observed in Ref. |2] can 
be associated to a conserved effective particle number for the low momentum 
regime. The strong turbulence solution predicts a characteristic dependence 
on the dimensionality of space. Our classical-statistical lattice simulations 
provide a striking confirmation of this dependence for three and four dimen- 
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Figure 2: The particle number spectrum of the d = 4 simulation for 
different times. The two graphs correspond to different initial conditions and 
different couplings as indicated. While the low-momentum scaling behavior 
is insensitive to these details, only the upper graph shows the expected 
Kolmogorov scaling at high momenta. 
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Figure 3: The particle number spectrum of the d 
O(10)-symmetric scalar field theory. 



3 simulation for the 
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sions. In particular, we see no indications for a dependence on the number 
of field components A^. This strongly suggests that the universal behavior 
associated to the conserved effective particle number indeed only depends on 
the dimensionality of space and the value of the dynamic scaling exponent 
z = 1 for relativistic dynamics with zero anomalous dimension. 
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